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Abstract 

In this paper we propose and analyze a new coupling procedure for the Hybridiz- 
able Discontinuous Galerkin Method with Galerkin Boundary Element Methods 
based on a double layer potential representation of the exterior component of the 
solution of a transmission problem. We show a discrete uniform coercivity estimate 
for the non-symmetric bilinear form and prove optimal convergence estimates for 
all the variables, as well as superconvergence for some of the discrete fields. Some 
numerical experiments support the theoretical findings. 
AMS Subject classification. 65N30, 65N38, 65N12, 65N15 

1 Introduction 

In this paper we propose and analyze a new coupling procedure for the Hybridizable 
Discontinuous Galerkin Method (HDG) [H] and a Galerkin Boundary Element procedure 
based on a double layer potential representation. The model problem is a transmis- 
sion problem in free space, coupling a linear diffusion equation with variable diffusivity 
(div kVm = /) in a polygonal domain, with an exterior Laplace equation. The transmis- 
sion conditions are given by imposing the value of the difference of the unknown and its 
flux on the interface between the interior and the exterior domains. 

Much has been written on the advantages and disadvantages of the many Discontin- 
uous Galerkin schemes. In support of this piece of work, let us briefly hint at some of 
the features that make the HDG -originally a derivation of the Locally Discontinuous 
Galerkin (LDG) method by Bernardo Cockburn and his collaborators- a family of inter- 
est. First of all, HDG uses polynomials of the same degree k > for discretization of all 
the variables -in the case of diffusion problems, a scalar unknown u, the vector valued flux 
q = —kVu and a scalar unkown on the skeleton of the triangulation-, attaining optimal 
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order 0{h^^^) in the approximation of all of them. This feature compares well with Mixed 
Finite Element Methods, of which HDG is a natural modification: HDG can be considered 
as a variant of the Raviart-Thomas and Brezzi-Douglas-Marini Mixed Elements, imple- 
mented with Lagrange multipliers on the interelement faces [T], eliminating the degrees 
of freedom used for stabilization and using a stabilization (not penalization) parameter 
to perform the same task. Similar to the Lagrange multiplier implementation of mixed 
methods by Arnold and Brezzi [1], HDG is implemented by hybridization, reducing the 
unknowns to those living on interelement faces. This reduces considerably the number 
of degrees of freedom and makes the method competitive with respect to other families 
of Finite Element Methods, while still being a mixed method, that approximates several 
fields simultaneously. Another interesting feature of HDG is the fact that for polynomial 
degrees > 1, the two unknowns related to the scalar field (the respective approximations 
inside the element and on the skeleton) superconverge at rate 0{h^^'^). This allows for the 
application of standard postprocessing techniques [7j that can be traced back to Stenberg 

Among many of its good properties (for the limited set of equations where it is us- 
able), the Boundary Element Method provides a reliable form of constructing high order 
absorbing boundary conditions with complete flexibility on the geometric structure of the 
domain (it does not require the domain to be convex or even connected). It is therefore 
natural to test and study the possibility of using BEM as a way of generating a coupled 
discretization scheme with HDG for transmission problems. The coupling of DG and 
BEM started with the study of LDG-BEM schemes [H], extended to other DG methods 
of the Interior Penalty (IP) family [llj. All of these methods, plus several new ones, 
were presented as particular cases of a methodology for creating symmetric couplings of 
DG and BEM in [10]. One of the methods in this latter paper, using HDG and BEM, 
was recently analyzed in [8]. As opposed to symmetric couplings with BEM, that need 
two integral equations (and thus four integral operators) to reach a stable formulation, 
non-symmetric couplings require only one integral equation (and two integral operators). 
From this point of view, they provide simpler (and more natural) forms of coupling BEM 
with field methods (FEM, Mixed FEM or DG). The non-symmetric coupling of Mixed 
FEM and BEM is recent [T7j and expands ideas used for analyzing the simple coupling 
of FEM-BEM [21]. The presentation and testing of non-symmetric coupling of DG meth- 
ods of the IP family with BEM [19], led to its analysis [15], providing the first rigorous 
proof of convergence of this kind of schemes. Note that, unlike in the case of symmetric 
couplings, where the stability analysis boils down to an energy argument, coercivity prop- 
erties in non-symemtric couplings pose a serious analytic challenge. The present paper 
makes a contribution in that direction, showing how to develop a coercivity analysis for 
non-symmetric coupling of HDG-BEM, and pointing out some interesting facts about the 
size of the diffusion parameter in the interior domain. 

The paper is structured as follows. In Section [2] we present the method. Sections |3] to 
[5] cover the analysis in progressive steps. Section |6] discusses several related extensions and 
possible modifications of the method. In Section[7]we include some numerical experiments. 
A collection of known results concerning transmission problems, potentials and integral 
operators, is given in Appendix |A] for ease of reference. 
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Regarding notation. Basic theory of Sobolev spaces is assumed throughout. Norms 
of L'^{0) type will be subscripted with the integration domain || ■ Ho. All other norms will 
be indexed with the name of the space. 



2 Formulation and discretization 

Let f2 C M*^ be a polygonal domain when = 2 or a Lipschitz polyhedral domain when d = 
3, with boundary F. For simplicity we will assume that F is connected. (Connectedness 
plays a minor role in the amount of energy-free solutions of the transmission problem.) 
Let fi+ := M'^ \ f2 be then the domain exterior to F. The outward pointing unit normal 
vector field on F is denoted n. 



2.1 Statement of the problem 



We consider a scalar positive diffusion coefficient n G L°^{Q) such that k ^ G 



(fi). 



We are interested in the following transmission problem: an elliptic second order diffusion 
equation in f2, written first order system, 



q + kVm = 



and 



div q = / in Q, 



(2.1a) 



the Laplace equation in the exterior domain 

Am+ = in and u = o(l) at infinity, (2.1b) 

coupled through two transmission conditions on the common boundary of these domains 



u 



u. 



+ /3o 



and 



q ■ n 



on F. 



(2.1c) 



In principle, boundary values of u and u+ in (2.1c) are taken in the sense of traces for 
functions with local Sobolev regularity. Similarly q • n and d-aUj^ are defined in a 
weak form as elements of the space H~^^'^{V). Basic regularity requirements for data is: 
/ G L2((]), /5o G H^'^{T), and /3i G ^^(F). We note that /3i G H-^'^{T) is valid data for 
the transmission problem, but not for the type of formulation we are going to use. Data 
functions are assumed to satisfy the following compatibility condition: 



/+ / /3 



0. 



(2.2) 



A more detailed discussion about this condition and solvability issues for the transmission 



problems is given in Section |A.4[ Here we just note the following: 



(a) When d 



2, condition (2.2) is necessary and sufficient for the transmission problem 



(2.1) to have a solution. In this case m+ = 0{r ) as r = |x| — )■ oo. 



(b) When = 3, problem (2.1) has a unique solution even if (2.2) fails to hold. 
(2.2) is satisfied, the solution decays as m+ = C(r~^) at infinity. Proposition A. 5 b) 



If 



shows the elementary modification of data that is needed to have condition (2.2) 
satisfied, while still providing a solution of the original problem. 
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2.2 Coupling with a double layer potential 

Condition ( |2.2[ ) allows us to write double layer potential 



(2.3) 



(See (A. 2) and (A. 6) for the definitions and Proposition A. 2 for the representation result.) 
The exterior trace and normal derivative of the double layer potential can be written in 
terms of boundary integral operators, as explained in Section A.3[ These operators are 
introduced in (A.12)-(A.14). They allow us to write the transmission conditions (2.1c) in 



the equivalent form 

u = \ip + ICif + and q ■ n = Wy? — Pi on V. 



Note that 



= and W<^ = q ■ n + /3i on F 



(2.4) 



(2.5) 



imply 



u{^,(t)) := (Wvp,0)r+ U> / = (q-n + /3; 



e H'^^{T). (2.6) 



At this moment, it is convenient to collect all the equations that make up the formulation 
we are going to discretize: 



(2.7a) 
(2.7b) 
(2.7c) 
(2.7d) 



The exterior field is reconstructed as a double layer potential m+ = Vip. The next lemma 





= 


in Q, 


divq 


= / 


in Q, 


u - {^(f + ICif) 




on r. 


(q ■ n, 0)r + a;(v?, 0) 


= {PiAW 


V0 e H^^\T) 



ensures that the integral boundary condition and the normalization condition (2.5) for 



the density (p are adequately encoded in the system (2.7). 



Lemma 2.1. If {q,u,(p) is a solution to (2.7) and the compatibility condition (2.2) is 
satisfied, then Jj.ip = 0. 



Proof. Integrating (2.7b) over Q and testing (2.7d) with = 1, it follows that 



/ /=(q-n,l)r = a;(^,l)-(/5i,l)r = |r| [ ^ - [ (3^, 
Jn Jr Jr 



because of Proposition |A.3[ This proves the result 



□ 
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2.3 Discretization with HDG and Galerkin BEM 



We start the discretization process by describing the discrete geometric elements. The 
domain Q is divided into triangles [d = 2) or tetrahedra [d = 3) with the usual conditions 
for conforming finite element meshes. A general element will be denoted K, Th will be 
the set of all elements and h := maxx^Th the maximum diameter of elements of the 
triangulation. The triangulation is assumed to be shape-regular. Extension of the forth- 
coming results to non-conforming meshes is relatively simple while not straightforward, 
requiring the use of some finely tuned results about HDG on general grids [3], 

We assume the existence of a triangulation of the boundary F, formed by line 
segments {d = 2) or triangles {d = 3). We write hr := maXegPh h^. From the point of 
view of analysis, it is immaterial whether this grid is related to Th or not. (Note that 
for the experiments we will only use meshes where is the trace of 7^ on F, which 
makes the implementation considerably simpler.) At the time of stating and proving the 
convergence theorems we will add some technical restrictions relating the mesh-sizes of 
the two grids. 

The set of edges {d = 2) or faces {d = 3) of all elements of the triangulation will be 
denoted £h and dTh = Ueg^-^e will denote the skeleton of the triangulation. Volumetric 
integration will be denoted with parentheses 

{u,v)k := / uv, (p,q);^:= / p • q, (u,v)r^ := y2{u,v)K, 
Jk Jk k^j-^ 

while integration on boundaries will be denoted with angled brackets 



{u,v)oK :-- 



dK 



dK\r- 



The unit normal outward pointing vector field on dK will be denoted ngx, or simply n, 
when there is no doubt on what the element is. 

Local discrete spaces will be composed of polynomials. The set Vk{K) contains all 
(i-variate polynomials of degree less than or equal to k and 7^k{K) := (Vk{K)Y. li e E £h 
or e G Th, 'Pfc(e) is the set of polynomials of degree less than or equal to k defined on e, 
i.e., the space of {d — l)-variate polynomials on local tangential coordinates. Four global 
spaces will be used for discretization: 



Vh 


:={v 






Vir G Th], 


Wh 


■.={w 




w\k e Vk{K) 


VA^ G Th}, 


Mh 


:={v: 


dTh K 


: v\e e Vk{e) 


Ve G £h}. 


Yh 


:={0 


F ^ M : 


0GC(F), 0|e 


eVk+iie) VeGF4 



(2.8a) 
(2.8b) 
(2.8c) 
(2.8d) 



The discretization method consists of using the Hybridizable Discontinuous Galerkin 



method (see [6], [7]) for equations (2.7a)-(2.7b)-(2.7c) (this part of the system can be un- 
derstood as an interior Dirichlet problem) and conforming Galerkin (Boundary Element) 



Method for equation (2.7d) (this part is considered as a hypersingular integral equation 
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for an exterior Neumann problem). We thus look for q/i G 14, G Wh, Uh G Mh, 
ifh e Yh satisfying 



and 



«^"^qh, r)r^ - (uh, div r)rft + {uh, r ■ n)aTu 


= 


Vr G Vh, 


(2.9a) 


-(q/,, Vw)ri, + (qh • n, w)ot,. 




G l^;,, 


(2.9b) 




= 


V?; G M,,, 


(2.9c) 


{uh, v)r - {\<^h + K^'^h, v)r 


= {l3o,v)r 


V?; G M,,, 


(2.9d) 


-{^h ■ n,0)r + w(v5/i,0) 


= (/3i,0)r 


V0 G 


(2.9e) 



Here qh is defined on the boundaries of the elements using the expression 

qh ■= qh\K + TaxiuhlK - Uh)naK ondK, 



(2.9f) 



where tqk '■ OK — t- M is a non-negative stabilization function that is constant on each 
edge/face of dK and such that it is strictly positive on at least one edge/face of each 
triangle. The stabilization function r can be double-valued on interlement edges/faces 
and q/j is in principle multiple valued, but its normal component is made to be single 



valued through equation (2.9c). Note also that equations (2.9c)-(2.9d) are tested with 
the same space, but because of the integration domain, they produce together as many 
equations as the dimension of M^. Using the local solvers for the HDG method [B], the 
system (2.9) can be reduced to a system with {uh,'^h) G Mh x Yh as the only unknowns. 



A discrete counterpart of Lemma 2.1 holds 



Lemma 2.2. If {qh,Uh,Uh,y^h) solves (2.9) and the compatibility condition (2.2) holds, 
then f^iph = 0. 



Proof. Testing equation (2.9b) with w 
the result follows. 



1 and equation (2.9e) and using Proposition A. 3 



□ 



The next sections deal with the analysis of the method: we first show the basic discrete 
coercivity arguments (Section [sj, proceed to proving energy estimates based on the HDG 
projection (Section |4]) and finally use duality arguments to analyze convergence (and 
superconvergence) of some of the fields (Section Isl). 



3 Solvability of the discrete system 



The aim of this section is to show that, under some conditions on k, the system (2.9) has 
a unique solution for h small enough. A relevant hypothesis will be the following: 

/^max := ||'^||L°°(n) < 4. (3.1) 

This hypothesis is related to coercivity of the underlying mixed formulation. We will 
discuss this hypothesis (and compare it with similar bounds for other coupled BEM-FEM 
schemes) in Section [7| where we will also make some experiments related to it. For some 
of the forthcoming arugments, we will use the piecewise constant function f) : dTh — ^ M 
given by f)|e := he- We will pay attention to the quantities 

Tmax := max ||r 11^00(9^) =: ||r||Loo(a-7;) and ||f)T||Loo(9rh)5 
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noting that because of the shape regularity of the grid, we can bound hK\\T\\L°°(^0K) < 
C||f)r 1 1^00(^-7^ ) for all K. Three more polynomial spaces will appear in our arguments: 



{u e Vk{K) : {u,w)k = G Vk-i{K)], 

{Vi{K)Y = {q G Vk{K) : (q,r)x = Vr G Vk-i{K)], 

{q:dK : qle Vk{e) We G S{K)}, 



where in the last space, we have denoted S{K) := {e & £h '■ e C dK}. 
Lemma 3.1. (a) If q E V^{K) and q = on e E £{K), then q = 0. 
(b) The following decomposition is orthogonal in L'^{dK): 

UkidK) = {^r\9K ■ n : V G Vi{K)} © {qU : q G Vi{K)}. 
Proof. Part (a) is straightforward. Part (b) is Lemma 4.1 of j9]. □ 
Lemma 3.2. Let Pk '■ L^(fi) — )■ Wh be the orthogonal projection onto Wh- Then 

\\u - PkuWdK < Ch][^\\Vu\\K WK G Th,Wu G H\n). 

Proof. Using the Bramble-Hilbert lemma (or a generalized Poincare inequality) , the trace 
theorem and an argument about finite dimensions, it is easy to prove that 

\\u-Pku\\QK<C\\Vu\\ji WueH\K), 

where K is the reference element and Pk '■ L'^{K) — )■ Vk{K) is the orthogonal projector. 
The result follows then by a scaling argument. □ 

Lemma 3.3 (A discrete coercivity estimate). For discrete triples (q^, Uh, Uh) G \4 x Wh x 

Mh satisfying 

{(liv(\h,w)T„ + {r{uh - Uh),w)aT,^ =0 Vw G Wh, (3.2a) 
(q/i ■ n + T{uh - Uh),v)aTh\T = WvE Mh, (3.2b) 

and G H^{Vt), we consider the quadratic form: 
Q ■= {K'^qh,(lh)rh + {T{uh-Uh),Uh-Uh)dTh + llVn^ll^ - {u^,(\h ■n + T{uh-Uh))v- 



If (3.1) holds, then there exist C > and D > such that if ||l)T||L°°(a7^,) ^ -D, 

Q > C(||q/i||n + llVu^lln + {T{Uh - Uh),Uh - Uh)dT^- 



Proof. Condition (3.2b) is equivalent to the fact that the discrete normal fluxes q/j • n + 
T{uh — Uh) are single valued on internal faces/edges. Therefore 

Eh ■■= {u^, q/i ■ n + T{uh - Uh))r = Eqk, where Eqk := (m*, q^ ■ n + T{uh - Uh))dK- 

K&Th 



If Pk '■ L^(fi) — )■ Wh is the orthogonal projection onto Wh (as in Lemma 3.2), then 

EdK = (qft ■ n, u^)aK + {r{uh - Uh), PkUi)aK + {T{uh - Uh),u^ - PkU^^dK 
= (q/i ■ n, u^)aK - (divq/,, PkU^)K + {'r{uh - Uh), - PkU^)dK 

= (q/i, Vu^)i^ + {T{Uh - Uh), - PkU^)dK, 



where we have used (3.2a) and the fact that divq/i G Vk-i{K). Therefore, by Lemma 3.2 
and Young's inequahty 



1^, 



dK\ 



< HqftlU'llVM^IlE- + C{hK\\T\\L^OK)f'^{T{Uh - Uh),Uh - M/,)^'!? || Vm^ || 



+ C'^||f)i"IU°°(arh)^ ^{r{uh - Uh),Uh - Uh)dK + 1{S + 
for arbitrary 6,6 > 0. Therefore 



K 



Q 



(k ^qh,qh) + {r{uh -Uh),Uh -Uh)dTh + II 



2 



E, 



h 



> (k 



> 



-1 
max 

1 - 

4 - K 



r^)||q.||E 



[1-C^\H)t\ 



K^ + ^))llv^.||^ 



L^ian)^ ^){T{uh - Uh),Uh - UhidTh 
A - K 



^max(4 ~l~ /^max 

where we have chosen 

r- 4 -|- Kmax 

This finishes the proof. 



|q/i|l5 + ^(^(w/i - Uh),Uh - Uh)dTh + 



6 :-- 



< 



16 



4 — 



max 1 1 T 1 1 9 

IIVm^L, 



16C2 



□ 



Proposition 3.4. // (3.1) holds, then there exists D > such that for all \\i)T\\L°°{dTh) — 
D, the system (2.9) is uniquely solvable. In particular, z/xmax < C , the system is uniquely 
solvable for h small enough. 

Proof. We only need to prove that the homogeneous system only admits the trivial solu- 
tion. Let then {qh,Uh,Uh,(Ph) be a solution of (2.9) with zero right hand side. By Lemma 



2.2 it follows that u{(ph,ij) = {'W(ph,4')r- We now test equations (2.9) with q^^, Uh, Uh, 



— q/i ■ n — T{uh — Uh) and Lph respectively, add the equations, and simplify the result to 
obtain 

(/t'^q,,, qh)Th + {'^{'^h-Uh),Uh-Uh)dTh + {-\Vh+^Vh, q/i-n+r(M,,-M/,))r+(>Vv2h, ^h)v = 0. 



Let now -u^ := T)(^h- Using then (A. 13) and Proposition A. 3, it follows that 
{n'^qh, qh)rh + {T{uh - Uh),Uh - Uh)dr,, + (7^^^*, q/i ■ n + r(u/, - Uh))r + || Vn* 



0, 



where 7~m* makes reference to the trace of uj\^ on F. We are now in the hypotheses of 
Lemma 3.3, which implies that qh = 0, Vm^ = 0, and T{uh — Uh) = on dK for all K . 
By Proposition A.l and the fact that Jp iph = 0, this implies that iph = 0. 



Going back to ( |2.9a ), it follows that 



{Vuh,r)K + {uh-Uh,r-n)9K = Mr eVk{K) WK 



(3.3) 
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and therefore 



By Lemma 3.1 ^b), there exists q G V^{K) such that Uh—Uh = q of OK. Since T{uh — Uh) 



and r 7^ on at least one e G £{K), it follows by Lemma [3A| ^a) that q = and therefore 



Uh = Uh on dTh- Using this information in (3.3) and testing with r = Vuh, it follows that 



Uh is piecewise constant. At the same time, Uh = Uh on dTh, which implies that both Uh 



and Uh are constant. Testing (2.9d) with v = 1 implies that this constant value has to 



vanish. This finishes the proof of uniqueness. □ 

4 Estimates by energy arguments 

Convergence analysis of the coupled HDG-BEM scheme follows the main lines of the 
projection-based analysis of HDG methods [7] (see also [S]). We start by recalling the 
local HDG projection [7]. Given (q, m), we define (Ilq, Hm) G T4 x Wh as the solution of 
the local problems 

(nq, r)^ = (q, r)^ Vr G Vk-iiK) (4.1a) 

(Uu, w)k = (m, w) Ww G Vk-iiK) (4. lb) 

{Hq- n + Tllu,v)dK = {q- n + TU,v)oK We7lk{dK), (4.1c) 

for all K E Th- We also consider the orthogonal projections P : riee^-fe -^^(^) ~^ and 
Pp : H^/\T) Yh. The error analysis is carried out by comparing the discrete solution 
with the projection of the exact solution, i.e., in terms of the quantities: 

el := Uq - qh, el := Uu - Uh, el := Pu - Uh, e^ := Pr^ - (fh- 



The error in the normal fiux (recall (2.9f)) is 

e/i := £^ ■ n + T{el — el) = n — qh ■ n, where vr := Ilq • n + t{IIu — Pu). 
Lemma 4.1. The following inequalities hold for all K and e G £{K): 

h,\\9-q-n\\l < c(||q-nq||i. + /i^||V(Pfcq-q)l|^), (4.2) 

||^-q-n||^ < ChK\\Vq\\l, (4.3) 

where P^q is the best L^(^7)'^ approximation of q in Vh- 

Proof. Note also that by definition of the projections 

{Ti,v)dK = {^-^,v)aK WeUkidK) \/K. (4.4) 

This shows that 7r|e is the best L'^{e) approximation of q ■ ng on Fk{e). Using then a local 
trace inequality, we can easily bound for e G S{K) 

he\\n - q ■ n\\l < h,\\Pkq - q\\l < c(\\Pkq - q\\l + hl\\V{Pkq - qjWl) ■ (4.5) 

The second inequality follows from a similar argument, comparing with Pgq instead of 
Pfcq. □ 
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For simplicity, in some of the forthcoming arguments we will shorten fC := \X + IC. 
The following three discrete functionals will be relevant in the sequel as well: 



ci(r) 

C2{y) 



i(nq-q),r)r;^, 
-(tt - q • n, 0)r + a;(Pr<^ - 



(4.6a) 

(4.6b) 
(4.6c) 



Proposition 4.2 (Energy inequality). // (3.1) holds, then there exists Hq such that for 
all h < Hq, 



(4.7) 



Proof. By the definition of the projections and (4.4), it follows that 



(nu,divr)7^^ + {Pu,r- n)or,^ 


= ci(r) 




Vr e Vh, 


(4.8a) 


-(nq, Vw)r^ + {n,w)ar^ 


= if, w)n 






(4.8b) 


-(7^>^)ar,\r 


= 






(4.8c) 


{Pu,v)r - {iCPr^,v)r 


= {1^0, v)r - 


C2(y) 


\fv G Mh, 


(4.8d) 


-(7f,0)r + w(PrV5,0) 


= (/3i,0)r + 


oM 


V0 G 


(4.8e) 



Subtracting (2.9) from equations (4.8) we obtain the error equations: 



-V„r)r,-(eLdivr)r,+ 


{el,r- n)ar. 


= ci(r) 


Vr G Vh, 


(4.9a) 




-£^),w)or. 


= 


G Wh, 


(4.9b) 




-{^h,v)dTh\'r 


= 




(4.9c) 




-{)Celv)r 


= -C2(y) 


yveMj,, 


(4.9d) 


-(4,0) 


r + w(e^,0) 


= C3(0) 


V0 G F;.. 


(4.9e) 



Testing these equations with e^, e'^, e^, —e^ and respectively, adding them, and sim- 
plifying, we obtain 



(4.10) 



Ci(s^J+C2(4)+C3(5D 



Let now e\ '.— I^Ef^ and note that Proposition |A.3 and ( |A. 13[ ) imply that 

2 



||V£^||rod\r + 



el 



and 



- = -7-4- 



Applying now the coercivity estimate of Lemma 3.3, it follows from (4.10) that 

+ {risl - el), el - ^l)^, + u:{elel) < c{c,{el) + c^ie^) + c,{el)). 



Since the bilinear form uj is -coercive (Proposition |A.3[), the result follows. □ 
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At this moment we have to introduce a technical hypothesis, 

hr<Cmm{hK:Kerh,Knr^(I)}, i.e., \\i)~^\\L^(r) < Ch^\ 



(4.11) 



meaning that the meshsize of Th is controlled by the minimum local meshsize of Th on 
r. This hypothesis is satisfied, for instance, if Th is quasiuniform in a neighbordhood of 
r and Th is equal the restriction of Th to F, or is a refinement of this inherited mesh. We 
will also assume that 



hr\\r\ 



L°°{r) 



< C. 



Theorem 4.3. Assume that (3.1), (4.11) and (4.12) hold. Then for ||[)T||Loc(g7-^) 
enough 



(4.12) 
small 



\4\\n + Hel - e%el - el)lff^^ + ||e^||^v.(r) < c(App^, + ||PrV^ 



where 



App^ := ||nq-q||n+ ($^4||V(Pfcq-q)||^ 

K 



1/2 



Proof. We will derive the proof from Proposition AJ2_ and some estimates related to the 
functionals defined in (4.6). First of all 

|ci(sDI<C||nq-q|b||s'Jb. 



(4.13) 



Using Proposition A. 4 and a scaling argument (note the use of (4.11) in the scaling 
argument), it follows that 



|c2(r,)| < ||ir(Pr¥^-^)||r(||4-n||r + +||r 

'1/2., gn II ||l/2 



1/2 



re. 



< C\\Pr^ - ^h[h-''''\\elU + Mll^rMel 



An Aubin-Nitsche duality argument and properties of best approximation by piecewise 
polynomial functions shows that 



\\P^ip - v^llr < Ch\!'^\\PYip - V5||i^i/2 



(r)- 



(4.15) 



For the third functional, we use Proposition A. 3 and (4.4) 

\cz{el)\ < \ - - n,el - Pel)r\ + C\\PrV - ^\\HU2^r)\\et\\m/2ir), (4.16) 

where it is to be understood that the local projection Pe:^ is only applied (and needed) 
on F. A local scaling argument (cf [131 Theorem 3.2]) shows that 



\r"\el-Pel)\\r<C\\el\\^. 



/2(r)- 



This inequality and Lemma 4.1 can then be used in (4.16) to prove that 



\cz{el)\ < C{k^^l + ||Pr<^ - <^||//i/2(r)) |k^||i/i/2(r) 



(4.17) 



(4.18) 



The result is then a consequence of Proposition 4.2 and the bounds (4.13), (4.14), (4.15), 
and (l4l8|. □ 
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Corollary 4.4. In the hypotheses of Theorem \4.!% 

||[)i/2(q ■ n - q, • n)|Ur. < C(App^ + ||Pr^ - ^||Hi/.(r))- 
Proof. Adding and subtracting vr, we can bound 

||[)i/2(q.n-q,-n)||ar. < ||[)^/V,||ar, + ||[)'/'(^ - q " n)||ar. 



A scaling argument (using the fact that G V^), plus Lemma 4.1 and Theorem 4.3 finish 
the proof. □ 

Let Fi, . . . , be the faces (with d = 3) or edges (when d = 2) oi T. We consider the 
space X"^(r) := nil ^""(r), endowed with its product norm. 

Corollary 4.5. Assume that the hypotheses of Theorem \4 .31 hold and that the exact so- 
lution of ^ satisfies: q G H^+^iVtY, u G H^+^iVt), ip G X^+^{T). Then 

||q-qftlb + y - '^h\\m/2(^r) 

< C/i^+"^(|q|^/fc+i(f^) + 6r\u\Hk+i[n)) + C/;,p^^^^||v2||x'=+2(r), 
where 6r '■= max;^ ||T||Loo(o_ft'\eK) ^^'^ ^ ^{^) such that Tf,^ = \\T\\Loo(^g^y 

Proof. This bound is a direct consequence of Theorem |4.3[ using well-known estimates 
for the best approximation by piecewise polynomials and [7^ Theorem 2.1]. □ 



5 Estimates by duality arguments 

Consider the problem 



K-^d + ve = 


in Q, 


(5.1a) 


div d = el 


in Q, 


(5.1b) 


Auj = 


in 


(5.1c) 


$(r) + C(r-'^+i) 


at infinity, 


(5. Id) 


e = uj 


on r , 


(5.1e) 


d • n = —dnUJ 


on r. 


(5.1f) 



By Proposition A. 2 and (5.1c)-(5.1f), we can represent u = VQ+S{d-n) in Therefore, 
using the definition of the operator W -see (A. 12)- and Proposition A. 4, we can write 

- (d ■ n, + /CV^)r + (WV^, e)r = G if^/='(r). (5.2) 

We will assume the following regularity estimate 

||0||/f2(Q) + ||d||Hi(n) < C'rcgllf^hlb- (5.3) 

Since this is a transmission problem, this hypothesis is mainly related to the regularity of 
the diffusion coefficient k. If k is smooth and equal to one in a neighborhood of F, then 
(5.3) holds (see [3 Proposition 3.4] for a similar argument). 
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Proposition 5.1. In the hypotheses of Theore'm \4.3\ and assuming the regularity estimate 
(5.3), it holds 

\Kh < Ch'^''^^'''^[Appl + \\Pr^ - ^||Hi/.(r)). (5.4) 



Proof. The gist of the proof consists of testing the error equations (4.9) with lid, — IIB, 
— PG , — (nd + r(nB — -PO)), and — Pr© respectively, add them and manipulate the 
result. Let us first test the first error equation (4.9a) with r = lid: 

{^-\q - q,), nd)r, - {el div nd)r, + (el, Hd ■ n)9r, = 0. 



Using integration by parts twice, the first condition in the definition of the HDG projection 

\\El\\l + {El{Ud-d).n)en. 



(4.1a) and the fact that divd = e^, we can write 
Adding these two equations, it follows that 



I u||2 



;q - q,), nd)7^^ - {el - el, (Hd - d) ■ n)ar, + (^1, d ■ n)&r. 



{K-\q - q,), nd)r, + {r{el - el), U& - P0)ar, + {el d ■ n)r, (5.5) 



where in the last equation we have applied (4.1c) and the fact that el and d are single 
valued on interelement faces. We next test (4.9b) with w = HQ and apply (4.1b) and the 
fact that Ve = -^-M: 

= {dwelQ)r, + {T{el~el),UQ)9r, 



t^-'el,d 



n,PQ)sr, + {T{el-el),UQ)9r,- 



The following step consists of subtracting this equation from (5.5) to obtain 



'hWn 



-i(q - q,), nd)r, - {^-'el d)r, + {^l d ■ n)r - (£^ n + r(e« - T^), Pe)ar. • 



(1) 



c 



(2) 



Bound for C'^\ Since n ^d = — V0, 



(5.6) 



(1) 



-\q - q^), nd - d)r, + (Hq - q, -n M)r, 

-i(q - q,), nd - d)r, + (nq - q, VG - Pfe-iV0)r„ 



where Pfc-i is the local orthogonal projection on the spaces 'Pk-i{K) if A; > 1 and P_i = 0. 
Note that the inclusion of Pfc_iVO is possible because of the definition of the HDG 
projection (4.1). Therefore 

\C'il\ < Ch\\q - q,|b(|d|Hi(n) + nu^m) + C/i--^'='i>||nq - q|b|eU2(^) 



< (||q - qftib + linq - q|b). 



(5.7) 



by regularity (5.3) and [71 Theorem 2.2]. 
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Bound for C^^ . Let now be the Clement approximation [5j| on a Pi conforming 
finite element space on a triangulation Th whose restriction to the boundary is Th- We 
then use (4.1c) to write 

{el, Tid)r - (q ■ n - q/j ■ n, Qh)r + (q ■ n - q/, ■ n, 6^ - P9)r, 



C, 



(2) 



(5i 



where 7?^ := lid ■ n + r(n6 — PQ). Testing the error equation (4.9d) with vf^, and 
using (5.2), it follows that 



{{\X+1C){^ - ^h).T^d - d ■ n)r + (>V(<^ - <^/,), e)r 
{{\X + /C)(<^ - <^,), TT, - d ■ n)r + (W(<^ - Vh), e - e;,)r 
+ (q-n-q/,-n,e,,)r, 



(5.9) 



where we have used that 



(q ■ n - q/, ■ n, 0)r + (>V(v3 - ^h), 0)r = ^{v - Vh, ' 



en. 



which is just (4.9e) after using Lemmas |2.1 and 2.2 Inserting (5.9) in (5.8) and applying 
Propositions A. 3 and A. 4 we can bound 

('2') ^ 

I < C\\(f-iph\\m/2(^r){\\TTd-d-n\\H~i/2(^j^) + \\Q-(S>h\\m/^r)) 



[)V2(q . n - q, . n)||r(||r^/'(i^e - e)||r + ||r^/'(e - e,)||r). (5.10) 



,-1/2 



-1/2/ 



We now just need to bound all the terms in the right hand side of (5.10). A duality 
argument (Aubin-Nitsche trick), Lemma 4.1, and the regularity assumption (5.3) show 
that 



l^d - d ■ n||j^-i/2(r) = ||P(d-n 

< Ch\\Vd\\n<Ch\\el\\n 



d • n|| j:^-i/2(p-) < Ch-p 



1/2, 



TTd - d ■ n| 



(5.11) 



Well-known properties of the Clement interpolant (see in particular [8, Proposition 5.2]) 
and (5.3) prove also that 

|e - e,||^i/2.n + llr'/'(0 - 0h)||r < Ch\e\H2ra) < Ch\\6l\\n. (5.12) 



Finally, we just overestimate 
||r'/'(P0-0)||r < ||r'/'(^fe0-0)l|r < C/i"''"^^''^|0| 



< Ch'^"'^^'''^\el\\n, (5.13) 



after applying a discrete trace inquality and (5.3). Bringing the bounds (5.11), (5.12) and 
(5.13) to (5.10) and using Corollary 4.4, we obtain the bound 



^(2)| < ^^min{l,.}||^«||^(App^ ^ _ + llp^y, _ 



(5.14) 



The bound (5.4) is now a direct consequence of (5.6), (5.7), (5.14) and Theorem 4.3 



□ 
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Note that in a bsence of any kind of additional regularity hypothesis, the estimate 
of Proposition 5.1 can be easily repeated without the additional /i™n{fc,i} ^j^^^^ provides 



superconvergence to the method. Also, as a consequence of Proposition 5.1, the local 
postprocessing technique of [H Section 5] can be applied here providing a Oih^^"^) ap- 
proximation of u. 



Corollary 5.2. In the hypotheses of Proposition [JJj for k >1, 



(r) 



(5.15) 



Proof. The argument to derive (5.15) from (5.4) can be taken verbatim from the proof of 
[Zt Theorem 4.1]. □ 



Proposition 5.3. In the hypotheses of Theorem \4.3\ 



II 



Proof. The error equation (4.9d) is equivalent to writing {}C{(p — iph),v)r = {e^,v)r for 



all V G Mh- Therefore, by Lemmas 2A and |2.2 and Proposition A. 4 , we can bound 

W^-^hh < C sup —- 

o^</,eHO(r) \m\T 

< C SUP {^k^lM^f^^c SUP ^^'-^'^^^ 
o^0eHO(r) ||0||r Q^4>eH0(r) \m\r 



< C 1 1 I j;^i/2(p) sup 

o^</.eHO(r) 



P0||j^-i/2(r) 



An Aubin-Nitsche duality argument shows then that ||0 — P4>\\H-^/^{r) 

< h^/^(p\\r- The 

proof then follows by Theorem 4.3 □ 



Corollary 5.4. In the hypotheses of Proposition 5.1 and assuming that h < Ch^, for 
k>l, 



Proof. It is a direct consequence of Proposition 5.3, using Corollary 5.2, the estimate 
^1"'/' and (1411]). □ 



l^llr ^ ||[) II L°°{T)\ 



6 Some additional considerations 

On the potential representation. The main drawback of a boundary integral for- 
mulation based on a potential ansatz on a polyhedral boundary is the expected lack of 
regularity of the associated density. This makes that the hypotheses for regularity of the 
solution in Corollary |4.5| might not be realistic. Note, however, that the numerical exper- 
iments shown in Section [7} for which the density is unknown, show that the computation 
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of the interior unknowns and of the exterior solution are not affected for the foreseeable 
lack of regularity of the density in the corners of the domain. On the other hand, the 
analysis allows for using a very refined grid Th near the corners. From the point of view 
of implementation, the case when Th is a refinement of the grid Th restricted to F, avoids 
many of the complications of dealing with general non-matching grids. 



Direct Boundary Integral formulation. A possible remedy for the above problem 
is the use of a direct formulation. For formulations of mixed type (and the one leading 
to the HDG method is one of such), this has been explained in great detail in [17]. The 



coupled formulation (2.7) has to be modified with the following arguments. First of all. 



we represent the exterior solution by 



where S is the single layer potential (|A.3|). Second, we consider the coupled problem 

(6.1a) 



K"^q + Vm° = 
div q = / 

u° - ^ = f3o 
~{q-n,l(f) + lC(t))T + w(v?, <P) = \^ 



+ lC<P)i 



in 
in 

on F, 

V0 G iJ^/2(F). 



(6.1b) 
(6.1c) 
(6. Id) 



Finally, the interior field u° is corrected by adding a constant, u = u° + where 



|r| 



((l,|V^ + /CV')r + (q-n + /3i,r/)r) 



7] 



$(| ■ -y|)dF(y), 



and $ is the fundamental solution for the Laplacian (A.l). This gives the solution to 
(2.1). Without the correction only q and u+ are correctly determined. If desired, it is 



possible to write ip + c = 7M+ and use this as a way of obtaining an approximation of 
the exterior trace. Equations (2.9) can be easily modified to handle this reformulation. 



From the point of view of implementation, this requires a simple rearrangement of the 



matrices in (2.9). The fact that the integral operators appear also in the right hand side, 
in a term of the form |0 + /C0)r, with G Y^^ can be easily handled by preprojecting 
the data function /3i in the space Mh restricted to the boundary F. Note that the HDG- 



BEM discretization of (6.1) leads to a system whose matrix is the transpose of the one in 



Section |2} Therefore, the quadratic form is the same, and all the energy arguments can be 
applied. The part of the analysis related to duality arguments requires some additional 
work though. 



Non-symmetric RT-BEM. A recent article |8j studies the symmetric coupling of 
HDG and BEM (using two integral equations) in parallel to the coupling of Raviart- 
Thomas mixed elements with BEM. The latter had appeared in the literature long ago 
[2]. Among other things, [8J provided an analysis of supercovergence in for the 
approximation of the variable m, and discussed the algorithmic advantages of hybridizing 
the RT method when coupled with BEM. Non-symmetric coupling of mixed elements with 
BEM has been proposed and studied in [17]. The analysis of this paper can be modified 
to include the methods proposed in [T7], thus providing some improved estimates that 
the variational techniques in that paper did not show. 
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On the diffusion parameter n. The final point for discussion is Hypothesis (3.1). 
First of all, let us mention that this inequality might not be sharp. We explore this in 
Section [7j However, there is evidence that some hypothesis like this is consubstantial to 
the mixed formulation that we are using for the HDG-BEM coupling. A similar result, 
with a lower bound for the minimum value of k -i.e., an upper bound for had 
been noticed in the context of FEM-BEM formulations: see [22j and [131 Theorem 5.1]. 
It has been recently proved jJOl Lemma 3.2] that there is actually a threshold for the 
minimum value of k under which the formulation used in non-symmetric FEM-BEM loses 
its well-posedness. Nevertheless, if the threshold is crossed in points at a certain distance 
of the coupling boundary, a compactness argument can be invoked to show that the 
formulation and its Galerkin discretizations (for sufficiently refined grids) are well posed. 
While global compactness arguments (where lower order terms are added and subtracted 
to the equation) do not seem to be applicable for the current format of HDG analysis, 
it is possible that some equivalent ideas could be used to prove that the methods of this 
paper can be used by placing the boundary F sufficiently far from the places where k is 
much larger than the exterior (unit) diffusivity. 



7 Experiments 

Convergence and superconvergence. For this example the domain is the square 
(0, 1) X (0, 1). The boundary mesh Th is taken to be the restriction of Th to F. The 
coarsest grid contains only two triangles. Other grids are obtained by uniform refinement, 
reaching up to 2048 elements and 3136 edges, of which 128 are boundary edges. Note that 
for a polynomial degree k, the dimension of the system that is solved (with Uh and (fh as 
unknowns) is -|- 1 times the number of edges plus the number of boundary edges. We 
take the diffusion parameter k{x, y) = 1 + x"^ and u{x, y) = exp(x + y) as exact solution 
in Q. The exterior solution is taken to be 

M+(x) := log y^^^ Xl = (0.3, 0.4), X2 = (0.7, 0.6). 

|X — X2 I 

The density (f is not known. To check errors we approximate the exterior solution at the 
observation point Xobs = (—0.1,0.1). Computation of the exterior solution is done using 
high order quadrature on the formula for the double layer potential. We tabulate and 
plot the following errors: 

q _ ||q - Qfelb u _ 11^ - UhWh + _ k+(Xobs) - M+.fe(Xobs)| 

l|q|b ' II^IU ' K(xobs)| 

||Hm-m/j|L \\Pu-Uh\\h 



e 



\u\\n \\u\\h 



The results are reported in Tables[T]to[3]for polynomials degree k = Oto k = 2 respectively. 
Estimated convergence errors are computed using consecutive grids. 
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ecr 




ecr 




ecr 


4 


ecr 




ecr 






Z.Oo 1 o{-l ) 




q 1 91 n/' 9^ 
0.1ZLU[-Z ) 








A 77/17/" 9^ 




2.2878(-l) 


0.77 


1.1908(-1) 


1.10 


1.4877(-2) 


1.07 


5.1969(-1) 


3.62 


2.6686(-2) 


0.84 


1.1867(-1) 


0.95 


5.5497(-2) 


1.10 


6.2955(-3) 


1.24 


3.5613(-1) 


0.55 


8.4902(-3) 


1.65 


6.0426(-2) 


0.97 


2.6790(-2) 


1.05 


3.5409(-3) 


0.83 


1.7215(-1) 


1.05 


3.8430(-3) 


1.14 


3.0435(-2) 


0.99 


1.3157(-2) 


1.03 


1.8998(-3) 


0.90 


8.6852(-2) 


0.99 


1.9418(-3) 


0.98 


1.5265(-2) 


1.00 


6.5192(-3) 


1.01 


9.8424(-4) 


0.95 


4.3598(-2) 


0.99 


9.9074(-4) 


0.97 



Table 1: Experiments for the lowest order {k = 0) method. All errors behave like 0{h). 





ecr 




ecr 




ecr 




ecr 




ecr 


9.1892(-2) 




3.2673(-2) 




4.8524(-3) 




9.3891(-2) 




9.9900(-3) 




2.5806(-2) 


1.83 


7.6059(-3) 


2.10 


7.6764(-4) 


2.66 


2.0486(-2) 


2.20 


1.6756(-3) 


2.58 


6.7443(-3) 


1.94 


1.7824(-3) 


2.09 


1.1487(-4) 


2.74 


1.2999(-3) 


3.98 


2.4053(-4) 


2.80 


1.7223(-3) 


1.97 


4.2887(-4) 


2.06 


1.5980(-5) 


2.85 


8.3227(-5) 


3.97 


3.2658(-5) 


2.88 


4.3549(-4) 


1.98 


1.0507(-4) 


2.03 


2.1137(-6) 


2.92 


1.1326(-5) 


2.88 


4.2930(-6) 


2.93 


1.0952(-4) 


1.99 


2.6000(-5) 


2.01 


2.7176(-7) 


2.96 


1.3269(-6) 


3.09 


5.5141(-7) 


2.96 



Table 2: Experiment for the case k = 1. Errors for q and u on dTh behave like 0{h'^). 
Comparison of u with respect to the projections (in elements and on their boundaries), 
as well as the exterior potential superconverge like 0{h^). 



Tests related to the diffusion parameter. For this example, the domain is the 
rectangle (—3/2,3/2) x (—1,1) and we use k = 0. We start with a fixed triangulation 
(produced with MATLAB's PDE Toolbox) with 936 elements. We take two different 
diffusion parameters: a constant value k(x) = Hcons and a piecewise constant function 



^int in (-3/4,3/4) x (-1/2,1/2), 
otherwise. 



We note that the values K.cons = and Hint = make the problem degenerate. We also note 
that the jump in the discontinuous diffusion coefficient is not resolved by the triangulation, 
i.e., we do not choose a triangulation with edges on the jump of the coefficient. As a test 
that might allow us to understand the effect of the diffusion parameter on the coupled 
system, we compute the condition number for increasing values of Hcons and and 
plot them together in Figure [TJ It is clear from the results for constant diffusion that the 



condition (4.1 ) is too restrictive -see also the comments in Section |6] concerning the recent 
results on nonsymmetric BEM-FEM-, but that growth of this parameter increases the 
condition number of the system significantly. It is also clear that if the diffusion coefficient 
grows far from the boundary, the problem is much better conditioned and that, for this 
case, the growth of conditioning appears to be linear in this parameter, which agrees with 
the basic fact that the matrix is an affine function of that seems not to degenerate as 
K grows. 

We then take several concrete values (three constant diffusion parameters and two 
piecewise constant, with the same notation as above), and plot the condition number for 
uniformly refined grids. Results are shown in Figure [2] 
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ecr 
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ecr 
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-4j 




1.7437(-3) 


2.85 


3.2768(-4) 


3.08 


3.9253(-5) 


3.92 


1.6677(-2) 


2.18 


9.6755( 


-5) 


3.28 


2.3006(-4) 


2.92 


3.8082(-5) 


3.11 


2.7636(-6) 


3.83 


5.8133(-5) 


8.16 


7.6652( 


-6) 


3.66 


2.8039(-5) 


3.04 


4.5361(-6) 


3.07 


1.4710(-7) 


4.23 


1.4569(-6) 


5.32 


4.1376( 


-7) 


4.21 


3.5088(-6) 


3.00 


5.5514(-7) 


3.03 


9.0340(-9) 


4.03 


2.6695(-9) 


9.09 


2.5795( 


-8) 


4.00 


4.3945(-7) 


3.00 


6.8673(-8) 


3.02 


5.7394(-10) 


3.98 


9.5947(-10) 


1.48 


1.6331( 


-9) 


3.98 



Table 3: Experiment for the case k = 2. Errors for q and u on dTh behave hke 0{h^). 
Comparison of u with respect to the projections (in elements and on their boundaries), 
as well as the exterior potential superconverge like The exterior field behaves 

somewhat erratically, which might be due to unaccounted errors in BEM quadrature. 



X 10 




Figure 1: Condition number of the system matrix for growing values of the diffusion 
parameter in the entire domain {kcous) and in an interior subdomain (Kint)- 



A Exterior and transmission problems: a compendium 

In this section we collect known results about layer potentials, as well as exterior and 
transmission problems associated to the Laplace equation. The results are straightforward 
consequences of results that are contained in [161 Chapters 6 & 8]. 

A.l Layer potentials 

Let 

^(r) - ^Jr) - I when d = 2, 
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Figure 2: Condition number of the system matrix for different diffusion parameters and 
uniformly refined grids. All of them behave like as could be expected from a two 

dimensional elliptic problem. 

be the fundamental solution of the Laplace equation. The double layer potential with 
density ip G H^/'^{T) is: 

(P^)(x) := ^Vy$(|x-y|)-n(y)(p(y)dr(y) 

= / 7 ' "i^^ ^(y) dr(y), vx G M'^ \ r. (a.2) 

The single layer potential is defined with a duality product: for A G i/^^/^(r), we define 
(5A)(x) := (A, $(|x - ■|))r, on any x G M'^ \ T. When A G L^(r), this duality product can 
be written in integral form 

(5A)(x) = ^$(|x-y|)A(y)dr(y). (A.3) 

Proposition A.l. Let (/? G H^^'^iV) and u := V(p. Then: 

u e H\n) and u e H\n+ n B{0; R)) Vi?, (A.4a) 

An = zn M'^ \ T, (A.4b) 

u G \ r) and u = C»(r-'^+^) as r = |x| ^ oo, (A.4c) 

— = (f and d^u — d^u = on T. (A.4d) 

Moreover, 

( —1 in / A ^\ 

Vl = -xn=< n n (A-5) 



0, m 

Proof. See [lEl Theorem 6.11] and [16. Chapter 8]. □ 
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The fact that constant densities produce a vanishing exterior solution of the Laplace 
equation motivates the introduction of the space 



hI'\T) := e H''\T) ■■1^ = 0}. (A.6) 

A. 2 Exterior solutions of the Laplace equation 

We consider functions u+ : — )■ M such that 

Am+ = in fi+ and u+ E H\n+ n fi(0; R)) Vi?. (A. 7) 

Since, by Weyl's lemma, locally integrable solutions of the Laplace equation are C°°, it then 
makes sense to impose a strong radiation condition at infinity. The general asymptotic 
condition we will deal with has the form 

u = Coo^{r) + 0{r^^'^) as r — )■ oo, uniformly in all directions. (A. 8) 

Note that this condition includes logarithmically growing solutions when d = 2, unless 
Coo = 0. The incoming flux on F is defined as 

Cfiux := -(<9n^i+, l)r- (A.9) 



Proposition A. 2. Let u satisfy (A. 7) and (A. 8). Then: 



(a) M+ admits the representation formula 

u+ = V-fu+ - SdnU+ = V{-fu+ + c) - SdnU+ Vc G Po(r). (A.IO) 

(b) Coo = Cflux and therefore, a necessary and sufficient condition for u to he decaying at 
infinity in the two dimensional case is caux = 0. 

1 /2 

(c) For every xq G fi, there exists a unique (p E Hq (T) such that 

M+ = Cflux$(|--xo|)+I^0. (A.ll) 

Therefore, can he represented as a douhle layer potential if and only if u+ = 
Q(j.-d+i^^ at infinity. 

Proof. Part (a) is a well-known representation formula Theorem 7.15]. The inclusion 



of any additive constant in the input of V follows from (A.5). Part (b) is straightforward 
using the strong integral form of the potentials. Part (c) can be easily proved by consid- 
ering u+ — Cflux '^'( I ■ — xqI) as the solution of an exterior Neumann problem and using a 
double layer potential representation [161 Theorem 8.19-8.21]. □ 
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A. 3 Integral operators 



Because of the transmission conditions satisfied by the double layer potential (A.4d), we 
can define the operators 

:= -d^V^, IC^ := U-f+V^ + I'V^). (A.12) 



Note that the conditions (A.4d) and the definition of /C imply that 

-f^Vip = ±lip + ICip. (A.13) 

When r is a polyhedral {d = 3) or polygonal {d = 2) boundary, the integral expression of 
the operator /C 

(/C¥^)(x) = \. / ^""7^^' Z^^^ V^(y) dr(y) (A.14) 



is valid on all points x G F that do not lie on edges. In particular (A.14) holds almost 
everywhere on F. Also, if F is a polyhedral boundary [d = 3) and (p,4' ^ L°°(F) are such 
that VrV^, Vr0 G L°°(F)'^ (here Vr is the tangential gradient), we have an integral form 
for the bilinear form associated to W: 

= / /(^W ^ Vrv^(x)) ■ (n(y) x Vr0(y)) $(|x - y |) dF(x)dF(y). (A.15) 
In the two dimensional case, the bilinear form is 



{Wif, <P) = dM^)dr(l>iy) $(|x - y|) dF(x)dF(y), (A.16) 
where dr is the tangential derivative on F. 

Proposition A. 3 (Properties of W). The operator W is bounded H^/'^{T) — )■ //"^/^(F). 
Its kernel is the set of constant functions. The bilinear form to : if^/^(F) x H^^'^(T) — t- M, 

1 /2 

is bounded, symmetric, and coercive. Also oj{ip, ip) = (Wf, f)r for all ip G Hq (F) and 

(Wip, ip)r = (Vm*, Vu^)^d\r, where = Dip, ip G iJ^/^(F). 

Proof See [ISl Theorems 8.20 & 8.21]. □ 

Proposition A. 4 (Properties of /C). The operator K, is bounded H^/'^iV) H^^'^iV) and 
H^iV) H^T). Moreover, 

{d^SX, ip)r = (A, -y + Mr VA G H-'/\T), ^ G H^'\V). 

Finally 

07^</)6L2(r) 

Proof. Boundedness in H^^'^iT) and the transposition property follow from the variational 
theory of layer potentials: see [HI Theorems 6.11 & 6.17]. For the if°(F) boundedness 
techniques of harmonic analysis are needed [21]. The bound (A. 17) follows from the fact 
that + /C is Fredholm of index zero and its kernel is the set of constant functions. □ 



U-^\U\\t<C sup ^^^"^.^^''^^^ Ve G H\T). (A.17) 
I'lJr o^</.6L2(r) ||0||r 
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A. 4 Transmission problems 

Let < K G L°°{Q) be such that G L°°{Q). The data of the transmission problem 
are / G L\n), /3o G //^/^(r), /3i G //-^/^(r). We look for m : ^ M and m+ : ^ M 



satisfying: the exterior Laplace equation (A. 7) with radiation condition (A. 8), the interior 
elhptic equation 

ueH\n) - V- (kVu) = / in fi, (A.18) 



and the transmission conditions 

■ju = 7M+ + /So on r 



and 



kVm ■ n = dnU+ + Pi on V. 



Because of (A.18)-(A.19), it is easy to prove that caux = Cdata, where 



Cdata '■= j f + l)r- 



(A.19) 



(A.20) 



Proposition A. 5 (Transmission problem). The transmission problem looking for u,u^ 



and Coo satisfying (A. 7), (A. 8), (A.18), and (A.19), is uniquely solvable 



(a) When d = 2, the solution is decaying (coo = in (A. 8) ) if and only if Cdata. = 0. 



(b) When d = 3, the solution satisfies u+ = 0{r ^) as r ^ oo if and only if Cd, 
Finally, if Cdata 7^ and xq G Q, we can write 

Cdata'^'d • -Xol) + U 



ata 



d 

+ ' 



where u'\. satisfies (A. 7), andu'\_ = 0{r ^) at infinity. The transmission conditions 
can then be written 



7M = 'ju'^ + /3o on r and 

where 

:= /3o + Cdata^d ■ -Xol), 

and therefore 



kVu ■ n = c^nwi + /3i on F, 



/3i:=/5i + CdataV$(|--xo|)-n, 
/ + (^i,l)r = 0. 



Proof. Existence and uniqueness of solution of (A. 7), (A. 8), (A.18), and (A.19) can be 



easily proved using a symmetric boundary- field formulation: see Section 1.5] for the 
general methodology and the needed background results. □ 



References 

[1] D. N. Arnold and F. Brezzi. Mixed and nonconforming finite element methods: 
implementation, postprocessing and error estimates. RAIRO Model. Math. Anal. 
Numer., 19(l):7-32, 1985. 



23 



[2] C. Carstensen and S. A. Funken. Coupling of mixed finite elements and boundary 
elements. IMA J. Numer. Anal, 20(3):461-480, 2000. 

[3] Y. Chen and B. Cockburn. Analysis of variable-degree HDG methods for convection- 
diffusion equations. Part I: General nonconforming meshes. IMA. J. Numer. Anal, 
(to appear). 

[4] Y. Chen and B. Cockburn. Analysis of variable-degree HDG methods for convection- 
diffusion equations. Part II: Semimatching nonconforming meshes. Math. Comp.. (to 
appear). 

[5] P. Clement. Approximation by finite element functions using local regularization. 
RAIRO Analyse Numerique, 9(R-2):77-84, 1975. 

[6] B. Cockburn, J. Gopalakrishnan, and R. Lazarov. Unified hybridization of discon- 
tinuous Galerkin, mixed, and continuous Galerkin methods for second order elliptic 
problems. SI AM J. Numer. Anal, 47(2): 1319-1365, 2009. 

[7] B. Cockburn, J. Gopalakrishnan, and F.-J. Sayas. A projection-based error analysis 
of HDG methods. Math. Comp., 79(271): 1351-1367, 2010. 

[8] B. Cockburn, J. Guzman, and F.-J. Sayas. Coupling of Raviart-Thomas and Hy- 
bridizable Discondinuous Galerkin Methods with BEM. SI AM J. Numer. Anal., 
50(5):2778-2801, 2012. 

[9] B. Cockburn and F.-J. Sayas. Divergence-conforming HDG methods for Stokes fiows. 
Math. Comp. (in revision). 

[10] B. Cockburn and F.-J. Sayas. The devising of symmetric couplings of boundary 
element and discontinuous Galerkin methods. IMA J. Numer. Anal., 32:765-794, 
2012. 

[11] G. N. Gatica, N. Heuer, and F.-J. Sayas. A direct coupling of local discontinuous 
Galerkin and boundary element methods. Math. Comp., 79(271):1369-1394, 2010. 

[12] G. N. Gatica and G. C. Hsiao. Boundary- field equation methods for a class of nonlin- 
ear problems, volume 331 of Pitman Research Notes in Mathematics Series. Longman, 
Harlow, 1995. 

[13] G. N. Gatica, G. C. Hsiao, and F.-J. Sayas. Relaxing the hypotheses of Bielak- 
MacCamy's BEM-FEM coupling. Numer. Math., 120(3) :465-487, 2012. 

[14] G. N. Gatica and F.-J. Sayas. A note on the local approximation properties of 
piecewise polynomials with applications to LDG methods. Complex Var. Elliptic 
Equ., 51 (2): 109-1 17, 2006. 

[15] N. Heuer and F.-J. Sayas. Analysis of a non-symmetric coupling of interior penalty 
DG and BEM. Math. Comp. (in revision). 



24 



[16] W. McLean. Strongly elliptic systems and boundary integral equations. Cambridge 
University Press, Cambridge, 2000. 

[17] S. Meddahi, F.-J. Sayas, and V. Selgas. Nonsymmetric coupling of BEM and mixed 
FEM on polyhedral interfaces. Math. Comp., 80(273):43-68, 2011. 

[18] S. Meddahi, J. Valdcs, O. Mcncndcz, and P. Perez. On the coupling of boundary 
integral and mixed finite element methods. J. Comput. Appl. Math.., 69(1):113-124, 
1996. 

[19] G. Of, G. J. Rodin, O. Steinbach, and M. Taus. Coupling of Discontinuous Galerkin 
Finite Element and Boundary Element Methods. SI AM J. Sci. Comput., 34-3:Al659- 
A1677, 2012. 

[20] G. Of and O. Steinbach. Is the one-equation coupling of finite and boundary element 
methods always stable? Berichte aus dem Institut fiir Numerische Mathematik, 
Bericht 2011/6, TU Graz, 2011. 

[21] F.-J. Sayas. The validity of Johnson-Nedelec's BEM-FEM coupling on polygonal 
interfaces. SI AM J. Numer. Anal, 47(5):3451-3463, 2009. 

[22] O. Steinbach. A note on the stable one-equation coupling of finite and boundary 
elements. SIAM J. Numer. Anal, 49(4):1521-1531, 2011. 

[23] R. Stenberg. Postprocessing schemes for some mixed finite elements. RAIRO Model 
Math. Anal Numer., 25(1):151-167, 1991. 

[24] G. Verchota. Layer potentials and regularity for the Dirichlet problem for Laplace's 
equation in Lipschitz domains. J. Funct. Anal, 59(3):572-611, 1984. 



25 



